---
title: "Replication code for 'Undeserving Heirs: How the Origins of Wealth Shape Attitudes Towards Redistribution'"
author: "Bastian Becker & Nora Waitkus"
date: "`r Sys.Date()`"
output: html_document
---


## Replication

The following code produces all results contained in 'Odious Wealth: How History Shapes Attitudes towards Redistribution'. Tables and figures are also saved to the folder 'output'.

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)

library(extrafont)
loadfonts(device = "win")
windowsFonts(Times = windowsFont("TT Times New Roman"))


library(lfe)
library(tidyr)
library(dplyr)
library(scales)
library(sjPlot)
library(ggplot2)
library(ggpubr)
library(RColorBrewer)
library(car)
library(texreg)
library(foreign)

load("data/data.Rdata")
```


```{r descriptives}

# Descriptive statistics ####

# balance plots
DAT <- DAT |>
  mutate(
    bal.right = (rile < 5) * 1,
    bal.left = (rile > 5) * 1,
    bal.eduhigh = edu_high,
    bal.old = (age > 52) * 1,
    bal.east = east,
    bal.female = female,
    bal.nonvoter = (party == "Non-voter") * 1,
    bal.spd = (party == "SPD") * 1,
    bal.cdu = (party == "CDU") * 1,
    bal.green = (party == "Green") * 1,
    bal.fdp = (party == "FDP") * 1,
    bal.afd = (party == "AfD") * 1,
    bal.theleft = (party == "Left") * 1,
    bal.other = (party == "Other party") * 1
  )

results <- DAT |>
  # Select only the computed variables and condition
  select(vgn_first, starts_with("bal.")) |>
  # Compute mean and standard error by condition for each variable
  group_by(vgn_first) |>
  summarise(across(
    starts_with("bal."),
    list(
      mean = ~ mean(.x, na.rm = TRUE),
      se = ~ sd(.x, na.rm = TRUE) / sqrt(sum(!is.na(.x)))
    ),
    .names = "{.col}_{.fn}"
  )) |>
  pivot_longer(-vgn_first,
               names_to = c("variable", ".value"),
               names_sep = "_")

results$variable <-
  factor(results$variable, levels = rev(unique(results$variable)))
levels(results$variable) <-
  rev(
    c(
      "Ideology: Right",
      "Ideology: Left",
      "Education: High",
      "Age: Old (53+)",
      "Place of living: East",
      "Gender: Female",
      "Vote Choice: Non-voter",
      "Vote Choice: SPD",
      "Vote Choice: CDU",
      "Vote Choice: Green",
      "Vote Choice: FDP",
      "Vote Choice: AfD",
      "Vote Choice: Left",
      "Vote Choice: Other"
    )
  )

ggplot(results, aes(x = variable, y = mean, col = vgn_first)) +
  geom_point(position = position_dodge(0.9)) +
  geom_errorbar(
    aes(ymin = mean - 1.96 * se, ymax = mean + 1.96 * se),
    width = 0.2,
    position = position_dodge(0.9)
  ) +
  labs(
    x = element_blank(),
    y = element_blank(),
    title = element_blank(),
    col = "Vignette"
  ) +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
        legend.position = "bottom") +
  coord_flip()

ggsave("output/app_balplot.pdf",
       width = 6,
       height = 4)

# Covariate balance tests

dat <- DAT |>  select(vgn_first, starts_with("bal.")) |>  na.omit()

mod0 <-
  lm(vgn_first == "Placebo" ~ 1, subset(dat, vgn_first %in% c("Placebo", "Nazi past")))
mod1 <-
  lm(
    vgn_first == "Placebo" ~ bal.other + bal.theleft + bal.afd + bal.fdp + bal.green +
      bal.cdu + bal.spd + bal.nonvoter + bal.female + bal.east + bal.old + bal.eduhigh +
      bal.left + bal.right,
    subset(dat, vgn_first %in% c("Placebo", "Nazi past"))
  )
anova(mod1, mod0)

mod0 <-
  lm(vgn_first == "Placebo" ~ 1, subset(dat, vgn_first %in% c("Placebo", "Apology")))
mod1 <-
  lm(
    vgn_first == "Placebo" ~ bal.other + bal.theleft + bal.afd + bal.fdp + bal.green +
      bal.cdu + bal.spd + bal.nonvoter + bal.female + bal.east + bal.old + bal.eduhigh +
      bal.left + bal.right,
    subset(dat, vgn_first %in% c("Placebo", "Apology"))
  )
anova(mod1, mod0)

mod0 <-
  lm(vgn_first == "Apology" ~ 1, subset(dat, vgn_first %in% c("Apology", "Nazi past")))
mod1 <-
  lm(
    vgn_first == "Apology" ~ bal.other + bal.theleft + bal.afd + bal.fdp + bal.green +
      bal.cdu + bal.spd + bal.nonvoter + bal.female + bal.east + bal.old + bal.eduhigh +
      bal.left + bal.right,
    subset(dat, vgn_first %in% c("Apology", "Nazi past"))
  )
anova(mod1, mod0)


# Outcome plots, all vignettes

dat <-
  as.data.frame(rbind(
    cbind(1, DAT$vgn_baseline_des),
    cbind(2, DAT$vgn_placebo_des),
    cbind(3, DAT$vgn_nazi_des),
    cbind(4, DAT$vgn_apology_des)
  ))
colnames(dat) <- c("Vignette", "Response")
dat$Vignette <- as.factor(dat$Vignette)
levels(dat$Vignette) <-
  c("Baseline", "Placebo", "Nazi past", "Apology")
dat$Response <- as.factor(dat$Response)
levels(dat$Response) <-
  c("Str Disagree", "Disagree", "Neither", "Agree", "Str agree")
dat <- aggregate(I(Response) ~ Vignette * Response, dat, length)
colnames(dat)[3] <- c("Frequency")


ggplot(dat, aes(y = Frequency, x = Vignette, fill = Response)) +
  geom_bar(position = "dodge", stat = "identity") + scale_fill_manual(values =
                                                                        c("grey80", "grey65", "grey50", "grey35", "grey20")) +
  xlab("") + theme_minimal()
ggsave("output/app_descriptives_des.pdf",
       width = 6,
       height = 2)

dat <-
  as.data.frame(rbind(
    cbind(1, DAT$vgn_baseline_red),
    cbind(2, DAT$vgn_placebo_red),
    cbind(3, DAT$vgn_nazi_red),
    cbind(4, DAT$vgn_apology_red)
  ))
colnames(dat) <- c("Vignette", "Response")
dat$Vignette <- as.factor(dat$Vignette)
levels(dat$Vignette) <-
  c("Baseline", "Placebo", "Nazi past", "Apology")
dat$Response <- as.factor(dat$Response)
levels(dat$Response) <-
  c("Str Disagree", "Disagree", "Neither", "Agree", "Str agree")
dat <- aggregate(I(Response) ~ Vignette * Response, dat, length)
colnames(dat)[3] <- c("Frequency")


ggplot(dat, aes(y = Frequency, x = Vignette, fill = Response)) +
  geom_bar(position = "dodge", stat = "identity") + scale_fill_manual(values =
                                                                        c("grey80", "grey65", "grey50", "grey35", "grey20")) +
  xlab("") + theme_minimal()
ggsave("output/app_descriptives_red.pdf",
       width = 6,
       height = 2)

# T-tests, baseline vs placebo vignette

t.test(DAT$vgn_placebo_des_num, DAT$vgn_baseline_des_num, paired = T)
t.test(DAT$vgn_placebo_red_num, DAT$vgn_baseline_red_num, paired = T)


# Outcome correlation, first vignette

cor(DAT$vgn_baseline_des_num, DAT$vgn_baseline_red_num)

```

```{r between}

# Between-subjects analysis ####

# Plot outcomes, first randomized vignette

DATb <-
  gather(DAT[, c("vgn_first", "vgn_first_des", "vgn_first_red")],
         condition,
         measurement,
         c(vgn_first_des, vgn_first_red),
         factor_key = TRUE)

DATb$condition <-
  ifelse(DATb$condition == "vgn_first_des",
         "Deservingness",
         "Redistribution")


DATb |>
  group_by(vgn_first, condition) |>
  summarise(mn = mean(measurement),
            se = sd(measurement) / sqrt(length(measurement))) |>
  mutate(vgn_first = factor(vgn_first, levels = c("Placebo", "Nazi past", "Apology"))) |>
  ggplot(aes(x = vgn_first, y = mn, fill = condition)) +
  geom_bar(position = position_dodge(), stat = "identity") +
  geom_errorbar(aes(ymin = mn - 1.96 * se,
                    ymax = mn + 1.96 * se),
                width = .2,
                position = position_dodge(.9)) +
  xlab(element_blank()) + ylab(element_blank()) +
  scale_y_continuous(limits = c(1, 5), oob = rescale_none) +
  scale_fill_manual(values = brewer.pal(4, "Paired")[c(4, 2)]) +
  labs(fill = "Outcome") +
  theme_minimal()

ggsave("output/main_betweensubjects.pdf",
       width = 6,
       height = 3)

# Pairwise t-tests, first randomized vignette

t.test(vgn_first_des ~ vgn_first, DAT[DAT$vgn_first %in% c("Placebo", "Nazi past"),])
t.test(vgn_first_red ~ vgn_first, DAT[DAT$vgn_first %in% c("Placebo", "Nazi past"),])
t.test(vgn_first_des ~ vgn_first, DAT[DAT$vgn_first %in% c("Apology", "Nazi past"),])
t.test(vgn_first_red ~ vgn_first, DAT[DAT$vgn_first %in% c("Apology", "Nazi past"),])

```

```{r within, results='asis'}

# Within-subjects analysis ####

# Regression analysis

DATw_des <-
  gather(
    DAT[, c(
      "caseid",
      "speeder",
      "right",
      "left",
      "edu_high",
      "age",
      "female",
      "east",
      "party",
      "vgn_placebo_des_num",
      "vgn_nazi_des_num",
      "vgn_apology_des_num"
    )],
    condition,
    measurement,
    c(vgn_placebo_des_num, vgn_nazi_des_num, vgn_apology_des_num),
    factor_key = TRUE
  )

DATw_des$condition <-
  recode(
    DATw_des$condition,
    "'vgn_placebo_des_num'='Placebo';'vgn_nazi_des_num'='Nazi past';'vgn_apology_des_num'='Apology'"
  )
DATw_des$condition <-
  factor(DATw_des$condition, levels = c("Placebo", "Nazi past", "Apology"))

mod1 <- felm(measurement ~ condition | caseid, DATw_des)
mod1bin <- felm(measurement > 3 ~ condition |
                  caseid, DATw_des)

DATw_red <-
  gather(
    DAT[, c(
      "caseid",
      "speeder",
      "right",
      "left",
      "edu_high",
      "age",
      "female",
      "east",
      "party",
      "vgn_placebo_red_num",
      "vgn_nazi_red_num",
      "vgn_apology_red_num"
    )],
    condition,
    measurement,
    c(vgn_placebo_red_num, vgn_nazi_red_num, vgn_apology_red_num),
    factor_key = TRUE
  )
DATw_red$condition <-
  recode(
    DATw_red$condition,
    "'vgn_placebo_red_num'='Placebo';'vgn_nazi_red_num'='Nazi past';'vgn_apology_red_num'='Apology'"
  )
DATw_red$condition <-
  factor(DATw_red$condition, levels = c("Placebo", "Nazi past", "Apology"))
mod2 <- felm(measurement ~ condition | caseid, DATw_red)
mod2bin <- felm(measurement > 3 ~ condition |
                  caseid, DATw_red)

# Regression tables

htmlreg(list(mod1, mod1bin, mod2, mod2bin))
  
texreg(
  list(mod1, mod1bin, mod2, mod2bin),
  digits = 3,
  sideways = F,
  caption.above = T,
  fontsize = "scriptsize",
  threeparttable = T,
  file = "output/app_results1.tex",
  booktabs = T,
  dcolumn = T,
  use.packages = F,
  caption = "Regression Results for Vignette Analysis",
  label = "table:results1",
  custom.gof.rows = list("Individual FE" = c("yes", "yes", "yes", "yes")),
  custom.header = list("Deservingness" = 1:2, "Redistribution" =
                         3:4),
  custom.model.names = c("(1) Cont. ", "(2) Dich.", "(3) Cont.", "(4) Dich."),
  custom.coef.map = list(
    "conditionNazi past" = "Nazi past",
    "conditionApology" = "Apology",
    "conditionPlacebo" = "Placebo"
  ),
  custom.note = "\\item\\textit{Note:} Linear regression (OLS) with individual fixed-effects; unit of analysis are vignettes (three per individual respondent); reference category = baseline vignette. Models 1 and 3 use continuous outcome variable, ranging from 1 (strongly disagree) to 5 (strongly agree); models 2 and 4 dichotomized outcome variables with 1 indicating (strong) agreement and 0 all other response options. (*=.05, **=.01, ***=.001)"
)

# Plot coefficients

plot_models(
  mod1,
  mod2,
  vline.color = "grey20",
  dot.size = 1.5,
  m.labels = c("Deservingness", "Redistribution"),
  legend.title = "Outcome",
  colors = brewer.pal(4, "Paired")[c(2, 4)]
) + theme(text = element_text(size = 20)) + theme_minimal()

ggsave("output/main_withinsubjects.pdf",
       width = 6,
       height = 2)

# Regression analysis, without speeders



mod1 <-
  felm(measurement ~ condition | caseid, subset(DATw_des, speeder == 0))
mod1bin <- felm(measurement > 3 ~ condition |
                  caseid, subset(DATw_des, speeder == 0))


mod2 <-
  felm(measurement ~ condition | caseid, subset(DATw_red, speeder == 0))
mod2bin <- felm(measurement > 3 ~ condition |
                  caseid, subset(DATw_red, speeder == 0))

# Tables

htmlreg(list(mod1, mod1bin, mod2, mod2bin))

texreg(
  list(mod1, mod1bin, mod2, mod2bin),
  digits = 3,
  sideways = F,
  caption.above = T,
  fontsize = "scriptsize",
  threeparttable = T,
  file = "output/app_results1_nospeeders.tex",
  booktabs = T,
  dcolumn = T,
  use.packages = F,
  caption = "Regression Results for Vignette Analysis (Excluding Speeders)",
  label = "table:results1_nospeeders",
  custom.gof.rows = list("Individual FE" = c("yes", "yes", "yes", "yes")),
  custom.header = list("Deservingness" = 1:2, "Redistribution" =
                         3:4),
  custom.model.names = c("(1) Cont. ", "(2) Dich.", "(3) Cont.", "(4) Dich."),
  custom.coef.map = list(
    "conditionNazi past" = "Nazi past",
    "conditionApology" = "Apology",
    "conditionPlacebo" = "Placebo"
  ),
  custom.note = "\\item\\textit{Note:} Linear regression (OLS) with individual fixed-effects, excluding the fastest 20\\% of respondents; unit of analysis are vignettes (three per individual respondent); reference category = placebo vignette. Models 1 and 3 use continuous outcome variable, ranging from 1 (strongly disagree) to 5 (strongly agree); models 2 and 4 dichotomized outcome variables with 1 indicating (strong) agreement and 0 all other response options. (*=.05, **=.01, ***=.001)"
)

# Within-subjects analysis, socioeconomic sub-groups ####

mod1_right <-
  felm(measurement ~ condition | caseid, subset(DATw_des, right==1))
mod1_left <-
  felm(measurement ~ condition | caseid, subset(DATw_des, left==1))
mod1_edulow <-
  felm(measurement ~ condition | caseid,
       subset(DATw_des, edu_high == 0))
mod1_eduhigh <-
  felm(measurement ~ condition | caseid,
       subset(DATw_des, edu_high == 1))
mod1_young <-
  felm(measurement ~ condition |
         caseid, subset(DATw_des, age <= 52))
mod1_old <-
  felm(measurement ~ condition | caseid, subset(DATw_des, age > 52))
mod1_west <-
  felm(measurement ~ condition |
         caseid, subset(DATw_des, east == 0))
mod1_east <-
  felm(measurement ~ condition |
         caseid, subset(DATw_des, east == 1))
mod1_men <-
  felm(measurement ~ condition | caseid,
       subset(DATw_des, female == 0))
mod1_women <-
  felm(measurement ~ condition | caseid,
       subset(DATw_des, female == 1))

htmlreg(list(
  mod1_right,
  mod1_left,
  mod1_west,
  mod1_east,
  mod1_edulow,
  mod1_eduhigh,
  mod1_young,
  mod1_old,
  mod1_men,
  mod1_women
))

texreg(
  list(
    mod1_right,
    mod1_left,
    mod1_west,
    mod1_east,
    mod1_edulow,
    mod1_eduhigh,
    mod1_young,
    mod1_old,
    mod1_men,
    mod1_women
  ),
  digits = 3,
  sideways = T,
  caption.above = T,
  fontsize = "scriptsize",
  threeparttable = T,
  file = "output/app_results1_subgroup.tex",
  booktabs = T,
  dcolumn = T,
  use.packages = F,
  caption = "Regression Results for Vignette Analysis: Deservingness by Subgroup",
  label = "table:results1_subgroup",
  custom.gof.rows = list("Individual FE" = rep("yes", 10)),
  custom.header = list(
    "Political Ideology" = 1:2,
    "Place of living" = 3:4,
    "Education" = 5:6,
    "Age" = 7:8,
    "Gender" = 9:10
  ),
  custom.model.names = c(
    "Right",
    "Left",
    "West",
    "East",
    "Low",
    "High",
    "Young",
    "Old",
    "Men",
    "Women"
  ),
  custom.coef.map = list(
    "conditionNazi past" = "Nazi past",
    "conditionApology" = "Apology",
    "conditionPlacebo" = "Placebo"
  ),
  custom.note = "\\item\\textit{Note:} Linear regression (OLS) with individual fixed-effects; unit of analysis are vignettes (three per individual respondent); reference category = placebo vignette. Deservingness as continuous outcome variable, ranging from 1 (strongly disagree) to 5 (strongly agree). (*=.05, **=.01, ***=.001)"
)

mod2_right <-
  felm(measurement ~ condition | caseid, subset(DATw_red, right == 1))
mod2_left <-
  felm(measurement ~ condition | caseid, subset(DATw_red, left == 1))
mod2_edulow <-
  felm(measurement ~ condition | caseid,
       subset(DATw_red, edu_high == 0))
mod2_eduhigh <-
  felm(measurement ~ condition | caseid,
       subset(DATw_red, edu_high == 1))
mod2_young <-
  felm(measurement ~ condition |
         caseid, subset(DATw_red, age <= 35))
mod2_old <-
  felm(measurement ~ condition | caseid, subset(DATw_red, age > 35))
mod2_west <-
  felm(measurement ~ condition |
         caseid, subset(DATw_red, east == 0))
mod2_east <-
  felm(measurement ~ condition |
         caseid, subset(DATw_red, east == 1))
mod2_men <-
  felm(measurement ~ condition | caseid,
       subset(DATw_red, female == 0))
mod2_women <-
  felm(measurement ~ condition | caseid,
       subset(DATw_red, female == 1))

htmlreg(list(
  mod2_right,
  mod2_left,
  mod2_west,
  mod2_east,
  mod2_edulow,
  mod2_eduhigh,
  mod2_young,
  mod2_old,
  mod2_men,
  mod2_women
))

texreg(
  list(
    mod2_right,
    mod2_left,
    mod2_west,
    mod2_east,
    mod2_edulow,
    mod2_eduhigh,
    mod2_young,
    mod2_old,
    mod2_men,
    mod2_women
  ),
  digits = 3,
  sideways = T,
  caption.above = T,
  fontsize = "scriptsize",
  threeparttable = T,
  file = "output/app_results2_subgroup.tex",
  booktabs = T,
  dcolumn = T,
  use.packages = F,
  caption = "Regression Results for Vignette Analysis: Redistribution by Subgroup",
  label = "table:results2_subgroup",
  custom.gof.rows = list("Individual FE" = rep("yes", 10)),
  custom.header = list(
    "Political Ideology" = 1:2,
    "Place of living" = 3:4,
    "Education" = 5:6,
    "Age" = 7:8,
    "Gender" = 9:10
  ),
  custom.model.names = c(
    "Right",
    "Left",
    "West",
    "East",
    "Low",
    "High",
    "Young",
    "Old",
    "Men",
    "Women"
  ),
  custom.coef.map = list(
    "conditionNazi past" = "Nazi past",
    "conditionApology" = "Apology",
    "conditionPlacebo" = "Placebo"
  ),
  custom.note = "\\item\\textit{Note:} Linear regression (OLS) with individual fixed-effects; unit of analysis are vignettes (three per individual respondent); reference category = placebo vignette. Redistribution as continuous outcome variable, ranging from 1 (strongly disagree) to 5 (strongly agree). (*=.05, **=.01, ***=.001)"
)

# Plot coefficients for sub-groups
set_theme()

p1 <- plot_models(
  mod1_right,
  mod1_left,
  mod2_right,
  mod2_left,
  vline.color = "grey20",
  dot.size = 1.5,
  axis.title = NA,
  m.labels = c("Right (Des.)", "Left (Des.)", "Right (Red.)", "Left (Red.)"),
  legend.title = NA,
  colors = "Paired"
) + ylim(c(-1, 1)) + guides(color = guide_legend(ncol = 2)) + theme_minimal() +
  theme(text = element_text(size = 20), legend.position = "bottom")


p2 <- plot_models(
  mod1_edulow,
  mod1_eduhigh,
  mod2_edulow,
  mod2_eduhigh,
  vline.color = "grey20",
  dot.size = 1.5,
  axis.title = NA,
  m.labels = c(
    "Non-Tertiary (Des.)",
    "Tertiary (Des.)",
    "Non-Tertiary (Red.)",
    "Tertiary (Red.)"
  ),
  legend.title = NA,
  colors = "Paired"
) + ylim(c(-1, 1)) + guides(color = guide_legend(ncol = 2)) + theme_minimal() +
  theme(text = element_text(size = 20), legend.position = "bottom")

p3 <- plot_models(
  mod1_young,
  mod1_old,
  mod2_young,
  mod2_old,
  vline.color = "grey20",
  dot.size = 1.5,
  axis.title = NA,
  m.labels = c("Young (Des.)", "Old (Des.)", "Young (Red.)", "Old (Red.)"),
  legend.title = NA,
  colors = "Paired"
) + ylim(c(-1, 1)) + guides(color = guide_legend(ncol = 2)) + theme_minimal() +
  theme(text = element_text(size = 20), legend.position = "bottom")

p4 <- plot_models(
  mod1_east,
  mod1_west,
  mod2_east,
  mod2_west,
  vline.color = "grey20",
  dot.size = 1.5,
  axis.title = NA,
  m.labels = c("East (Des.)", "West (Des.)", "East (Red.)", "West (Red.)"),
  legend.title = NA,
  colors = "Paired"
) + ylim(c(-1, 1)) + guides(color = guide_legend(ncol = 2)) + theme_minimal() +
  theme(text = element_text(size = 20), legend.position = "bottom")

p5 <- plot_models(
  mod1_men,
  mod1_women,
  mod2_men,
  mod2_women,
  vline.color = "grey20",
  dot.size = 1.5,
  axis.title = NA,
  m.labels = c("Men (Des.)", "Women (Des.)", "Men (Red.)", "Women (Red.)"),
  legend.title = NA,
  colors = "Paired"
) + ylim(c(-1, 1)) + guides(color = guide_legend(ncol = 2)) + theme_minimal() +
  theme(text = element_text(size = 20), legend.position = "bottom")


# Add titles
p1 <- p1 + labs(title = "(a) Political Ideology") +
  theme(plot.title = element_text(hjust = 0.5, size = 20))
p4 <- p4 + labs(title = "(b) Place of Living") +
  theme(plot.title = element_text(hjust = 0.5, size = 20))
p2 <- p2 + labs(title = "(c) Education") +
  theme(plot.title = element_text(hjust = 0.5, size = 20))
p3 <- p3 + labs(title = "(d) Age") +
  theme(plot.title = element_text(hjust = 0.5, size = 20))
p5 <- p5 + labs(title = "(e) Gender") +
  theme(plot.title = element_text(hjust = 0.5, size = 20))


### Arrange and export
zp1 <- ggarrange(p1, p4, p2, p3, p5, 
                 ncol = 2, nrow = 3,
                 hjust = rep(-1, 5), align = "v",
                 font.label = list(size = 20, face = "plain"))

# Export
ggsave("output/main_withinsubjects_combined.pdf",
       width=12,
       height=8)


# Within-subjects analysis, party sub-groups ####

mod1_party <- list()
mod2_party <- list()
pl <- list()
i = 0
ttl<-c("(a) SPD","(b) CDU/CSU","(c) Green Party","(d) FDP","(e) AfD", "(f) Left Party"," (g) Other Party", "(h) Non-voter")

for (p in c("SPD","CDU","Green","FDP","AfD","Left","Other party","Non-voter")) {
  i <- i + 1
  mod1 <-
    felm(measurement ~ condition | caseid,
         subset(DATw_des, party == p))
  mod2 <-
    felm(measurement ~ condition | caseid,
         subset(DATw_red, party == p))
  
  mod1_party[[i]] <- mod1
  mod2_party[[i]] <- mod2
  
  print(
  pl[[i]] <- plot_models(
    mod1,
    mod2,
    vline.color = "grey20",
    dot.size = 1.5,
    m.labels = c("Deservingness", "Redistribution"),
    legend.title = "Outcome",
    colors = brewer.pal(4, "Paired")[c(2, 4)]
  ) +
    scale_y_continuous(limits = c(-1.25, 1.25)) +
    theme_minimal() + ylab("") +
    theme(text = element_text(size = 20), plot.title = element_text(hjust = 0.5, size = 20)) + 
    labs(title = paste(ttl[i]) 
          )
  )
}


### Arrange and export
zp2 <- ggarrange(plotlist = pl, 
                 ncol = 2, nrow = 4,
                 hjust = rep(-1, 5), align = "v",
                 font.label = list(size = 20, face = "plain"), 
                 common.legend = T, legend = "bottom")

# Export
ggsave("output/main_withinsubjects_party_combined.pdf",
       width=12,
       height=8)

# Regression tables

htmlreg(mod1_party)

texreg(
  mod1_party,
  digits = 3,
  sideways = T,
  caption.above = T,
  fontsize = "scriptsize",
  threeparttable = T,
  file = "output/app_results1_party.tex",
  booktabs = T,
  dcolumn = T,
  use.packages = F,
  caption = "Regression Results for Vignette Analysis: Deservingness by Vote Choice",
  label = "table:results1_party",
  custom.gof.rows = list("Individual FE" = rep("yes", length(levels(
    DAT$party
  )))),
  custom.model.names = levels(DAT$party),
  custom.coef.map = list(
    "conditionNazi past" = "Nazi past",
    "conditionApology" = "Apology",
    "conditionPlacebo" = "Placebo"
  ),
  custom.note = "\\item\\textit{Note:} Linear regression (OLS) with individual fixed-effects; unit of analysis are vignettes (three per individual respondent); reference category = placebo vignette. Deservingness as continuous outcome variable, ranging from 1 (strongly disagree) to 5 (strongly agree). (*=.05, **=.01, ***=.001)"
)

htmlreg(mod2_party)

texreg(
  mod2_party,
  digits = 3,
  sideways = T,
  caption.above = T,
  fontsize = "scriptsize",
  threeparttable = T,
  file = "output/app_results2_party.tex",
  booktabs = T,
  dcolumn = T,
  use.packages = F,
  caption = "Regression Results for Vignette Analysis: Redistribution by Vote Choice",
  label = "table:results2_party",
  custom.gof.rows = list("Individual FE" = rep("yes", length(levels(
    DAT$party
  )))),
  custom.model.names = levels(DAT$party),
  custom.coef.map = list(
    "conditionNazi past" = "Nazi past",
    "conditionApology" = "Apology",
    "conditionPlacebo" = "Placebo"
  ),
  custom.note = "\\item\\textit{Note:} Linear regression (OLS) with individual fixed-effects; unit of analysis are vignettes (three per individual respondent); reference category = placebo vignette. Redistribution as continuous outcome variable, ranging from 1 (strongly disagree) to 5 (strongly agree). (*=.05, **=.01, ***=.001)"
)

# General attitudes ####

dat <-
  aggregate(I(reparation_max) ~ nazi_economy + reparation_max,
            DAT,
            sum)
colnames(dat)[3] <- c("Frequency")
dat$nazi_economy <- as.factor(dat$nazi_economy)
levels(dat$nazi_economy) <-
  c("0-99", "100-199", "200-299", "300-399", "400-500")
DAT2 <- dat
DAT2$out <- "Maximum"
colnames(DAT2)[2] <- "reparation"

plot1 <-
  ggplot(dat, aes(x = nazi_economy, y = reparation_max)) + geom_point(aes(size =
                                                                                    Frequency), color = "grey70") +
  geom_smooth(
    data = DAT,
    aes(x = nazi_economy, y = reparation_max),
    method = "lm",
    fill = NA,
    color = "black"
  ) +
  xlab("Perceived Nazi Heritage") + ylab("Reparative Measures (max)") +
  theme_minimal() + theme(axis.text = element_text(size = 7))


dat <-
  aggregate(I(reparation_n) ~ nazi_economy + reparation_n, DAT, sum)
colnames(dat)[3] <- c("Frequency")
dat$nazi_economy <- as.factor(dat$nazi_economy)
levels(dat$nazi_economy) <-
  c("0-99", "100-199", "200-299", "300-399", "400-500")


plot2 <-
  ggplot(dat, aes(x = nazi_economy, y = reparation_n)) + geom_point(aes(size =
                                                                                  Frequency), color = "grey70") +
  geom_smooth(
    data = DAT,
    aes(x = nazi_economy, y = reparation_n),
    method = "lm",
    fill = NA,
    color = "black"
  ) +
  xlab("Perceived Nazi Heritage") + ylab("Reparative Measures (n)") +
  theme_minimal() + theme(axis.text = element_text(size = 7))

ggarrange(
  plotlist = list(plot1, plot2),
  ncol = 2,
  common.legend = T,
  legend = "bottom"
)

ggsave("output/main_descriptives_collab.pdf",
       width = 6,
       height = 3)

# Regression analysis, Composite index

mod0 <-
  lm(nazi_economy ~ female + age + east + edu_high + rile, DAT)
mod1 <-
  lm(reparation_n ~ nazi_economy + female + age + east + edu_high + rile,
     DAT)
mod2 <-
  lm(reparation_max ~ nazi_economy + female + age + east + edu_high + rile,
     DAT)

htmlreg(list(mod0, mod1, mod2))

texreg(
  list(mod0, mod1, mod2),
  digits = 3,
  sideways = F,
  caption.above = T,
  fontsize = "scriptsize",
  threeparttable = T,
  file = "output/main_results2.tex",
  booktabs = T,
  dcolumn = T,
  use.packages = F,
  caption = "Regression Results for Perceived Extent of the Nazi Economic Heritage",
  label = "table:results2",
  custom.header = list("Nazi heritage" = 1, "Reparative Measures" = 2:3),
  custom.model.names = c("(1) ", "(2) Number", "(3) Max."),
  custom.coef.map = list(
    "nazi_economy" = "Nazi heritage",
    "female" = "Female",
    "age" = "Age",
    "east" = "East-Germany",
    "edu_high" = "University degree",
    "rile" = "Ideology"
  ),
  custom.note = "\\item\\textit{Note:} Linear regression (OLS) with individual observations. Nazi heritage indicates estimate of the 500 largest companies with a Nazi past, from 1 (0-99) to 5 (400-500) (*=.05, **=.01, ***=.001)"
)

# Regression analysis, Individual items


mod1 <-
  lm(
    reparation1 == "Yes" ~ nazi_economy + female + age + east + edu_high + rile,
    DAT
  )
mod2 <-
  lm(
    reparation2 == "Yes" ~ nazi_economy + female + age + east + edu_high + rile,
    DAT
  )
mod3 <-
  lm(
    reparation3 == "Yes" ~ nazi_economy + female + age + east + edu_high + rile,
    DAT
  )
mod4 <-
  lm(
    reparation4 == "Yes" ~ nazi_economy + female + age + east + edu_high + rile,
    DAT
  )
mod5 <-
  lm(
    reparation5 == "Yes" ~ nazi_economy + female + age + east + edu_high + rile,
    DAT
  )


htmlreg(list(mod1, mod2, mod3, mod4, mod5))
  
texreg(
  list(mod1, mod2, mod3, mod4, mod5),
  digits = 3,
  sideways = F,
  caption.above = T,
  fontsize = "scriptsize",
  threeparttable = T,
  file = "output/app_results3.tex",
  booktabs = T,
  dcolumn = T,
  use.packages = F,
  caption = "Regression Results for Perceived Extent of the Nazi Economic Heritage II",
  label = "table:results3",
  custom.header = list("Reparative Measures" = 1:5),
  custom.model.names = c(
    "(1) Apology",
    "(2) Memory",
    "(3) Voluntary reparation",
    "(4) Mandatory rep.",
    "(5) Expropriation"
  ),
  custom.coef.map = list(
    "nazi_economy" = "Nazi heritage",
    "female" = "Female",
    "age" = "Age",
    "east" = "East-Germany",
    "edu_high" = "University degree",
    "rile" = "Ideology"
  ),
  custom.note = "\\item\\textit{Note:} Linear regression (OLS) with individual observations. Nazi heritage indicates estimate of the 500 largest companies with a Nazi past, from 1 (0-99) to 5 (400-500) (*=.05, **=.01, ***=.001)"
)

```
